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EFFECT OF AN APERTURE ON MEASUREMENT OF THE AXIAL 
DISTRIBUTION FUNCTION IN A MAGNETICALLY 
CONFINED PLASMA 
by Roman Krawec 
Lewis Research Center 

SUMMARY 

Theoretical calculations have been performed to determine the distortion in the 
axial velocity distribution when a magnetically confined plasma passes through an 
aperture. A general solution is obtained in terms of aperture length to radius ratio, 
transverse to axial temperature ratio, and normalized length. Results are shown for 
magnetic field strengths ranging from zero to infinity. It is found that the use of strong 
magnetic fields allows the aperture to pass a plasma sample from which the true axial 
distribution function may be obtained. 

The effects of distortion in the axial velocity distribution on a measurement of axial 
temperature by taking the slope of the natural logarithm of current versus voltage is also 
discussed. It is shown that the error in such a measurement can generally be kept below 
30. 5 percent. 


INTRODUCTION 

The usual methods of finding the axial energy distribution of particles within a mag- 
netically confined plasma consist of allowing a small sample of the plasma to pass 
through an aperture and performing an energy analysis on the resulting low-density 
plasma. This energy analysis may be done by using electrostatic fields alone (refs. 1 
and 2) or by means of a combination of electric and magnetic fields (refs. 3 and 4). The 
aperture may be placed at the end of a mirror machine (refs. 1 and 4), may be attached 
to the body of a probe which is placed into the plasma (ref. 2), or can consist of a mag- 
netically shielded duct (ref. 3) so that the particles are extracted transverse to the 
magnetic field. 



Analytical solutions of the change in the distribution function of particles passing 
through the aperture have generally been restricted to treating the particles as having 
rectilinear motion (refs. 1 and 2), in which case the effects of the magnetic field are 
considered negligible. Another method, which also neglects magnetic fields, has been 
to consider the aperture as a boundary separating two regions of different electric field 
strength and then determining what effect these fields have on the distribution function 
(refs. 3 and 5). 

One exception is a recent paper by Anderson, Eggleton, and Keesing (ref. 6), who 
treat the case of a point source of plasma in a magnetic field and show that the particle 
distribution is strongly affected when the plasma passes through an infinitely thin aperture. 

It is the purpose of this report to extend the calculations of Anderson, Eggleton, 
and Keesing to the case of a uniformly distributed, magnetically confined plasma with 
anisotropic velocity distribution which is allowed to pass through apertures of arbitrary 
length and diameter. The general method of solution is followed by applications to 
apertures of specific length and diameter. These solutions are then used to propose 
criteria for aperture design which will give minimum distortion in the axial distribution 
function. 


THEORY 

Formulation of Problem 

Consider a flat plate of thickness L immersed in a uniform magnetic field which is 
normal to its surface. Let the region to one side of this plate contain a plasma of known 
velocity distribution and the other side be evacuated. The plate is considered to contain 
a cylindrical aperture of radius R through which a portion of the plasma may flow. The 
situation to be considered is depicted in figure 1. 

The primary goal of this report is to find the axial velocity distribution of the parti- 
cles emerging from the aperture while a secondary goal is to determine whether an 
accurate measurement of axial plasma temperature can be obtained at the aperture exit. 
The variables to be considered are the aperture dimensions, the magnetic field strength, 
and the axial and transverse particle temperatures. The following assumptions are made: 

(1) Any particle striking the bounding surfaces of the aperture is absorbed. 

(2) The plasma is collisionless; that is, the mean free path is large compared with 
the aperture dimensions. 

(3) The regions under consideration are free of electrostatic fields. 

(4) Particle absorption at the walls is the only loss mechanism. 
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Magnetic field 

Figure 1. - Schematic of situation to be analyzed. 


Because of the presence of a magnetic field, the charged particle trajectories are 
helices spiraling around the magnetic lines of force. It is therefore convenient to re- 
solve this helical motion into a circular motion about the particles' guiding center super- 
imposed on the motion of the guiding center along the direction of the magnetic field. 

This direction is taken as the z-axis. In the remainder of this report, the term particle 
orbit, or simply orbit, will denote only this circular motion projected on the entrance 
plane. 

While the final objective is to calculate the number of particles emerging from the 
aperture per unit time with axial velocities between v_ and v + dv , the analysis will 
first separate the particles into two classes depending on whether the gyroradius is less 
than or greater than R. A transmission function for the flow of particles through the 
aperture will be calculated as a function of gyroradius and axial velocity for each class 
and the results will then be combined to compute a distribution function for the particles 
emerging from the aperture in terms of the initial distribution assumed. 

Particles entering the plane of the aperture with gyroradii in the range 0 < r g <R 
are placed in class I; while particles with gyroradii in the range R < r^ < °° are placed 
in class n. Class I particles can further be divided depending on whether their orbits 
are wholly within (Class la) or partly within (Class lb) the aperture. The classes are 
summarized below. 


Class 

Range on r 

& 

Description 

la 

lb 

n 

0 ^ r < R 

0S rg<H 
Rs r g < °° 

Orbit completely inside of aperture 
Orbit partly intersects aperture 
Orbit partly intersects aperture 
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The procedure followed herein will be to determine a transmission function which 
can then be multiplied by a velocity distribution function and integrated to obtain the 
velocity distribution function at the aperture exit or which can be multiplied by the 
velocity distribution function times the velocity and integrated to determine a particle 
flow at the aperture exit. The same results could be obtained by initially treating this as 
a flow problem and calculating particle flow directly. The particular approach using a 
transmission function was considered more useful. 

Consider a uniform distribution of particles with a given gyroradius r . The loca- 
tion of the guiding centers of these particles on the x-y plane will also be uniformly dis- 
tributed. Consider the particles within a plane slice of thickness dz such that the sur- 
face density of their guiding centers is ct. Now select an area dA in the plane of the 
aperture such that any particle of gyroradius r whose guiding center lies within dA 

o 

will intersect the aperture with some or all of its orbit. The total number of such parti- 
cles will be adA. If the length of the particle orbit projecting across the aperture is S, 
the probability of that particle entering the aperture as it reaches the aperture entrance 
plane is S/27rrg. The quantity crSdA/27rrg. then represents the total number of particles 
with gyroradii r^ which enter the aperture. Since the quantity S varies throughout A, 

this expression must be integrated. The total number of particles of gyroradii r which 

/ » & 

S dA with limits of integration appropriate to the 

value of r being considered. This area will be either circular or annular. In terms of 

the radial distance from the center of the aperture p, the integral becomes 




Sp dp 


( 1 ) 


The limits of integration in this expression depend on the class to which the particles 
have been assigned. 

If the particles belong to class I (r < R) 

o 



The integral whose limits are 0 to R - r g corresponds to particles belonging to class la 
while the remaining integral corresponds to particles belonging to class lb. 
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If the particles belong to class IE (r ^ R) 

o 


/r g+ R 

Sp dp 
-R 


( 3 ) 


Calculation of a Transmission Function 

While traversing the aperture length L, each of the particles will complete a number 
of orbits which depend on the strength of the magnetic field, the aperture length, and the 
axial velocity of the particle. The orbital angle through which each particle rotates 
while traversing the aperture is given by 


0 = 


m v„ 


(4) 


where all the quantities used are defined in appendix A. 

Therefore, the particles will travel through an arc of length r g 0 during the time it 
takes to pass through the aperture. Referring to figure 2, not all the particles which 
enter the aperture along the arc S will pass through, but only those which enter the 
aperture on that portion of the path given by S - 0r g . Consequently, the equations cor- 
responding to equations (2) and (3) for the exit plane of the aperture will be given by 


N, 


thru = ff0 < R 


.. /’^lu 

O + j- j (S 

g */R-r„ 


Or) p dp 


(5) 


for class I (r < R) and 

o 



(S - 0r g )p dp 


( 6 ) 


for class II (r ^ R). The limits have been changed from those used in equations (2) 
and (3) to insure that the quantity S - dr will never be negative since negative values of 

D 
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Figure 2. - Schematic depicting the portion of 
its orbit (denoted by heavy line) that a 
particle must be on in order to pass through 
the aperture. 


the integrand represent those particles that will impinge on the aperture walls. These 
limits can be rigorously specified in terms of r^ and v z , as shown in appendix B. 

In the first class (r < R), the integrand would range from 277 at the lower limit 
to 0 at the upper limit if 6 were zero. For any given value of 9, the upper limit is 
effectively reduced. The integral vanishes if 9 exceeds 277. 

In the second class (r ^ R), the integrand is zero at both limits when 9 is zero. 
The effect of increasing 9 is to shrink both limits until at last they coalesce. 

Since for any value of r , the uniformity of spatial distribution implies that the total 

& 2 
number of particles which will arrive at the aperture entrance should be cnrR , a trans- 
mission function for the aperture can be defined by dividing the expressions given by 

o 

equations (5) and (6) by (ottR ) (see appendix C). This results in 



In order to carry the calculations any further, we need the relation between S, p, 

and r . Referring to figure 3, 

& 


S = <pr 
Y g 


( 8 ) 
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Figure 3. - Schematic representation of rela- 
tion between S, p, and r g . 


The law of cosines gives 


from which 




2r p 
g^ 


S = 



where 


a, = 

2r 

g 


a 2 " 


< R2 - 


2r 


g 


( 9 ) 


( 10 ) 


(Ha) 


(lib) 


Hence, if we momentarily neglect the limits, the problem of finding a closed form 
expression for the transmission function consists of evaluating an integral of the form 




- e 


P dp 


( 12 ) 
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As shown in appendix C, the closed form expression for this integral is given by 


p£=<p 2 cos ~ 1 i a ip- — 


a 2\ ep 2 + 1 + 4a i a 2 sin -i ( 2a i p2 ’ 1 " 2a l a 2 


4a‘ 




V 1 + 4a l a 2 

£-a 2 p 4 + (1 + 2a x a 2 )p 2 - a|j ^ ^ 


lb 


2a. 


( 13 ) 


Hence, the transmission function becomes 


'R - r 


V = 


R 


£ ) + -M i C u 

' ttR 2 R_r g 


(r*. < R) 


(14a) 


-L. [if 2 * 

77R 2 Pl 


(r g ~ R ) 


(14b) 


Using such a transmission function, it is then simple to take a known velocity dis- 
tribution, for example, a two -temperature Maxwellian, and determine its shape at the 
exit of the aperture. 


Velocity Distribution at Aperture Exit 

Let the particles in the bulk plasma have a velocity distribution function given by 


F i< v x’V v z) = 


n Q m 


m \ 1/2 m(v 2 +v 2 )/2kT i -mv 2 /2kT z 

1 e y e 


2jrkT 1 l2ffkT 


(15a) 


n„m / ™ \1/ 2 -mv 2 /2kT . -mv 2 /2kT 

F i< v l» v z) = — e v l e 

112 kT ± l2?rkT J 1 


(15b) 
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The axial velocity distribution ahead of the aperture is then given by 


F i< v z> - 


n 0 m /_m_\ 1/2 e - mv z/ 2kT z 

kT ± UjrkT J 


/' 


-mv,/2kT . 
e 1 1 v ± dv x 


„ / m \V2 -mv 2 /2kT z -mv 2 /2kT z 

= ° fer 1 


( 16 ) 


The corresponding expression for the exit of the aperture is 


„ -mvj?/2kT 

W = r^- A i e Z Z 


kT, 


/ 


-mvJ/2kT i 


v ± t? dv ± 


(17) 


Current at Aperture Exit 

The arrival rate of particles at the aperture entrance gives rise to a current which 
can be expressed as 


dl i (v z ) = ttR qAj e 


2 _ . -”V 2kT z 


v„ dv„ 
z z 


(18) 


At the aperture exit, the corresponding expression is 


d If(v z ) = 


ffR^qmAj^ -mv z /2kT,. 


kT, 


(/' 


-mvf/2kT . 

e 1 v^ dvj dv z (19) 


If it were not for the effect of the aperture, the current could be analyzed in terms of 

some retarding potential which cuts off particles having velocities less than v = 2qV/m. 

max 
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In this case equation (18) can be expressed as 


M v z > v max) = nR 


■/' 


-mv?/2kT 


Z Z J 

e v z dv z 


max 


2 / kT z\ V2 ^ V / kT z 
- " rV »( ^ • 


( 20 ) 


The maximum current density occurs when V = 0 and is just the term in front of the 
exponential in equation (20). 


V v z ^ v max^ ^max exp ^ kT ^ 


( 21 ) 


From which 


qV _ 

kT r 


= In W - In > v max» 



In Ii<v z 


v max^ 



(22a) 


(22b) 


Equation (22b) is the widely used relation which permits T_ to be determined from 
the slope of a retarding -potential curve. The question arises as to whether this same 
procedure is applicable to the aperture exit. Unfortunately, a closed -form integration 
is not possible for this case. One is required to resort to numerical methods to obtain 
a set of retarding potential curves for a range of values of R, L, T L , T z , and B. It 
appears possible to condense the requirements by a normalization process. Suppose we 
normalize velocities in terms of their average values, the aperture radius in terms of 
an average gyroradius, and aperture length in terms of the distance that a particle of 
average axial velocity will travel while completing one cyclotron orbit. This leads to 
the following normalized parameters: 
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rv 

*2kT_ z 
z 

(23a) 

rv 

V 1 =1 


(23b) 

r>v 

R - 

qBR 

(23c) 

XV — 

^/mkTj^ 

rv 

L = - 

qBL 

(23d) 


2ff|/2mkT z 


It was also found convenient to normalize 


S 


and p as follows: 


r^J 



_S 

R 


R 

This permits rewriting the equations in the form 



2jrR 2 kT 


m 


•/ 


F f (v z )v z dv 2 


max 


(24a) 


(24b) 


(25a) 


F f (v z ) = n Q 



rv rv rv 


?=(1 

B n 1-r 


S 


<? g < i) 


(25b) 


(25c) 
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= -[iL 2u 

n p 7 


(? g - 1) 


(25d) 


1 = 


Q -1 / P 2 + ? ? - 

p COS 


2 p r 


g 


9 / o 2 ~2 1 1 

-ge! < -sin- 1 ( p g 

2r 

g 


j 4 /?f + l' 


'2,2 


4 + 2 


■2 <»-? 


IP - 




1/2 


(25e) 


max 


m /^ 

¥ kT z 


(25f) 


The limits of integration on p can be expressed as 






P lu = r g cos -r- sin 


? sin 2 * 

s \ 2 > 


(26a) 


= r cos 


p 2u~ 



? 2 sin 2 m 


(26b) 


ro 


p l = "g C ° S U 


H ~2 • 2 6 
'1 - r„ sin / - 


(26c) 


The following relations are useful: 


L_ = _0_ 
v z 2 tt 


(27a) 


rv 

v t 


r_ = 


R 


(27b) 


The use of relations (27) allows the expressions for current and the distribution 
function to be set up for computer integration for a range of values of L, T /T and 
L/R. 
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Special Case B ■ 0 


At first glance, it might appear that the case B = 0 is treated by simply examining 
equations (25) to (27) in the limit as B goes to zero. This, however, can lead to 
erroneous conclusions because of the extensive manipulation involving quantities which 
either go to zero or become unbounded when the magnetic field is set equal to zero. 

We thus proceed to rederive the problem for the case when the magnetic field is 
absent. The similarities to the case when a magnetic field is present will be pointed 
out as they arise. 

A particle with axial velocity v_ will traverse the aperture length in a time t 

Zl 

given by 


r = 


L 


v 


z 


(28) 


During this time, the particle will also travel a transverse distance 6, given by 


6 = V ± T = 



(29) 


Note that the quantity 6 corresponds to the arc length r 9 for the case B * 0 and, 

D 

in fact, is identical to it in value. 

Suppose the position of the particle as it enters the aperture is taken as the center 
of a circle of radius 6. This circle defines the possible positions that every particle 
of velocities v and v. can have after traversing an axial distance L. (The initial 
position of the particle plays the role of the guiding center for the previously treated 
case, while the circle which defines the location of a particle at the aperture exit cor- 
responds to the particle orbit. ) Following the procedure established previously, we ask 
what portion of this circle of radius 6 lies within the aperture. The entering particles 
are again broken up into two major classes depending on the value of 6, which are 
summarized as follows: 


Class 

Range on 6 

Circle of radius 6 

la 

Jb 

n 

0 < 6 < R 
0 < 6 < R 
R < 6 ^ 2R 

Lies completely within the aperture 
Lies partly within the aperture 
Lies partly within the aperture 
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(a) Class Ia ; all particles enter. (b) Class Ib ; d2n enter. 



(c) Class II; o/2tt enter. 


Figure 4. - Representation of the fraction of particles in each class 
capable of leaving the aperture for the case B = 0. 

Referring to figure 4, the fractions of particles capable of leaving the aperture are 
given by 



A transmission function (similar to eq. (7)) can now be defined for the zero field case 
and is given by 


2 f -i (pq + S 2 - R 2 \ 

2p dp + — - / cos — p 0 dp 0 




r R + - r 2 \ 

L cos 
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where (0 < 6 < R) 



(30b) 


where (R < 6 < 2R). The integral in equation (30b) is similar to one previously used and 
can be expressed in closed form. This gives 


cos 


- U " I) + 




(31a) 


where (0s^< l). 



(31b) 


where (1 < £ < 2) and 



R Rv z 


Thus, if 0 ^ | ^ 2 



The axial distribution function at the exit plane may now be written as 


F 0 (v z> = n 0 



(32) 


(33) 


(34) 
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where 


v = 


2\/2Rv z 


L 



(35) 


The current -voltage characteristic is then found in the usual manner. 


Limit as Magnetic Field Approaches Infinity 

The limiting case of arbitrarily large magnetic fields may be treated by noting that 
the orbital angle 9 always becomes greater than its maximum permissible value of 2jr 
and the gyroradius goes to zero as the magnetic field approaches infinity. An examination 
of equations (25c) and (26a) reveals that the transmission function is identically equal to 
1, in which case 


F (v z ) = 

rear 



B— 0 


= n 0 




(36) 


This equation is immediately recognizable as the axial distribution function at the 
entrance plane of the aperture. This is as it should be since the gyroradii of all particles 
are zero, and therefore the particle motion is purely along the magnetic field lines. Thus 
every particle that enters the aperture will reach the exit plane. 


RESULTS AND DISCUSSION 

Equation (25b) was set up for machine integration in terms of the normalized radius 
R, the length to radius ratio L/R, and the transverse to axial temperature ratio T 1 /T z . 
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The parameters L/R and T./T were varied between values of 0. 1 and 10, and K 
was varied over a sufficient range so that the curves approached the values taken on at 
B = 0 and B — °°. The case B = 0 (eq. (34)) was also calculated, the parameter of interest 
being (L/R) ^T/T. 

The computer programs used to perform the calculations are included in appendix D. 
The distribution functions were normalized by dividing through by the factor 
nQym/27J’kT z , and some typical results are presented in figure 5. 

The effects of changing the normalized radius are presented in figure 5(a); this is 

no 

equivalent to changing the magnetic field since R is directly proportional to the strength 
of the magnetic field. Specifically 


R = 420 ° - 00 B R (electrons) 
ml/2 

1 L 


(37a) 


9 800 B R 

t 1/2 

L 1 


(protons) 


(37b) 


where B is given in tesla, R in meters, and Tj^ in electron volts. 

Looking further at figure 5(a) indicates that the distribution function lies near the 
B — oo case when R » 10 and lies near the B = 0 case when R « 1.0. In general, a 
large value of R means that the aperture is much larger than the radius of gyration of 
most of the particles (strong field case), and thus the distribution function will remain 
relatively undisturbed. A small value of R implies that the motion of the particles 
will be nearly rectilinear as they pass through the aperture (weak field case). In this 
case, the distribution function will take on the characteristics of the B = 0 case. 

Figures 5(b) and (c) give the results of changing Tj/T z and L/R. The rather 
complex effects obtained by varying either of these parameters are clearly indicated. 

Further insight into the effects of varying either the length to radius or the tempera- 
ture ratio may be gained by looking at some limiting cases of B = 0. These are pre- 
sented in figure 6 for various values of the combined parameter (L/R)|/T 1 /T z . For a 
given value of this combined parameter the distribution function for nonzero magnetic 
fields will fall between the curve for (L/R) |/T 1 /T z equal to zero and the curve cor- 
responding to the given value of this combined parameter. 

As previously mentioned, the distribution function for large values of magnetic field 
will be unchanged from that at the aperture entrance. Noting first the case 
(L/R) yT/T z = 0. 1, it is clear that this curve lies quite close to the undisturbed distri- 
bution function for most values of the dimensionless axial energy and departs from it 
only near the origin. (The case (L/R) |/T i /T z = 0 is the undisturbed distribution function. 
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(a) Effect of varying the magnetic field (normal- 
ized radius). Length to radius ratio, 2.0; 
temperature ratio, 1.0. 


(b) Effect of varying the transverse to axial temper- 
ature ratio. Length to radius ratio, 2.0; norma] 
ized aperture radius, 0.889. 



(c) Effect of changing the length to radius ratio. 
Temperature ratio, 1.0; normalized aperture 
radius, 0.889. 


Figure 5. - Typical distribution functions at the aperture exit showing the effects of varying the magnetic field, the tem- 
perature ratio, and the length to radius ratio. 
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Figure 6. - Change in normalized distribution 
function du e to va riation of combined para- 
meter L/R^T | /T Z for the case B = 0. 

This is also the limiting case as the magnetic field approaches infinity for all values of 
L/R and T./T since particles move parallel to the B -field in this case.) Looking at 
the other curves shows that increasing either L/R or T./T_ has two effects 

(1) The curve representing the distribution function moves further away from the 
undisturbed case. 

(2) The portion of the curve which is nonlinear on a semilogarithmic plot occupies a 
larger range of values of dimensionless energy. 

The most severe region of distortion is always near the origin and is a maximum for 
large values of T./T . All the distribution functions shown will reach a region of con- 
stant slope for sufficiently large values of the normalized axial energy, where a true 
temperature may be measured. This region may never be reached in practice. The 
maximum value of the normalized energy reached will depend not only on the initial 
particle current available but also on how quiescent the plasma in question is. Generally, 
it may not always be practical to measure a current versus voltage curve over a range of 
variation of current greater than 100. Taking these restrictions into account, tempera- 
tures were evaluated by drawing a best fit straight line through the logarithmic current 
versus voltage curve (on a plot of the distribution function such as fig. 6) in the region 
of normalized energy from 4 to 5 and using its slope to determine a temperature. Using 
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Figure 7. - Ef fect of the combined para- 
meter L/r/Tj /T z on measured axial 
temperature. "Magnetic field strength, 0. 


a range of values of normalized energy closer to zero will give rise to greater errors in 
temperature measurement, while moving the region further away from zero will give rise 
to smaller errors. 

Figure 7, shows the ratio of the temperature that will be indicated by a measurement 
to actual temperature for B = 0. The indications are that the measured temperatures 
will be too large unless L/R is restricted to values less than 0. 1. The measured tem- 
peratures depart more and more as L/R is increased, reaching a maximum value of 
1. 305 times the actual temperature. Continued increases in L/R have no further effect 
on this measured temperature. 

The variation of measured temperature with magnetic field (normalized radius) is 
the subject of figure 8. It is clear that increasing L/R has the effect of bringing the 
point where the measured temperature is equal to the true temperature nearer to the 

Length to 
radius 



Normalized radius, R 


Figure 8. - Effect of normalized radius on a measurement 
of axial temperature. (The parameter R is directly 
proportional to the strength of the magnetic field. ) 
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origin. This indicates that comparatively small values of magnetic field will be effective 
in obtaining an undistorted value of the temperature for large values of L/R. 

The preceeding discussion gives preliminary criteria for aperture design (i.e., either 
L/R must be small or the magnetic field must be chosen so that the aperture is operated 
in the strong field region). 

In either case, it is not considered good practice to make the radius of the entrance 
apertures larger than the Debye length of a particle within the plasma. To establish an 
upper limit to collector current, consider an aperture with radius equal to a Debye length. 
The collected current is 


nv 2 

I = q — ttR^ (38) 

4 

where v„„ is the average axial velocity. With R equal to the Debye length, aperture 

cLV 

current is found to be independent of density and is 

I = 6. 58 T juA (39) 

z 

for electrons when the axial temperature is given in electron volts. 

An aperture which is built following these suggested criteria should yield a plasma 
efflux that truely represents the distribution function. 


SUMMARY OF RESULTS 

Expected changes in the axial distribution of velocities have been calculated for a 
collisionless plasma immersed in a magnetic field when such a plasma is allowed to flow 
through an aperture of finite length and diameter. The velocity distribution at the aper- 
ture exit was obtained in terms of the length to radius ratio of the aperture, the trans - 
verse to axial temperature ratio, and the normalized aperture length. 

The distribution functions for selected values of the above parameters were pre- 
sented along with the weak and strong magnetic field limits. The ratio of aperture length 
to radius and the ratio of transverse to axial temperature were varied from 0. 1 to 10. 0, 
while the magnetic field was varied from zero to infinity. The following results were 
noted: 

1. For very strong magnetic fields, the distribution function remains undistorted re- 
gardless of the values of the other parameters. 
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2. M inim um distortion in the case of weak or intermediate values of magnetic field 
occurs when the aperture length to radius ratio is very much less than one. 

3. The maximum ratio of temperature indicated at the aperture exit to that at the 
aperture entrance was 1. 305. 

Lewis Research Center, 

National Aeronautics and Space Administration, 

Cleveland, Ohio, January 6, 1970, 

129-02. 
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APPENDIX A 


SYMBOLS 


Aj constant, defined in 

eq. (16) 

dA element of area 

a^ defined in eq. (11a) 

a ~2 defined in eq. (lib) 

B magnetic field strength 

Ff(v z ) normalized axial velocity 

distribution function at 
aperture exit 

Fj(v x , v , v z ) velocity distribution func- 

tion at aperture 
entrance 

F.(v ) axial velocity distribution 

function at aperture 
entrance 

F.(v, , vj. , velocity distribution func- 
i' l 7 z'front 

tion at aperture 
entrance 

F()(v z ) axial velocity distribution 

function at aperture 
exit (B = 0) 

I integral, defined in 

eq. (12) 

no 

I integral, defined in 

eq. (25e) 

Ij(v z ) current at aperture 

exit 

L(v z ) current at aperture 

entrance 


M*z > v max> 

current at aperture en- 
trance due to particles 
whose axial velocities 
are greater than v w> „„ 

III d A 

*max 

maximum current at 
aperture entrance 

k 

Boltzmann constant 

L 

aperture length 

no 

L 

normalized aperture 
length 

m 

particle mass 

N 

total number of particles 
of gyroradius r^ en- 
tering aperture 

n o 

particle number density 

q 

electronic charge 

R 

aperture radius 

R 

normalized aperture 
radius 

r s 

gyroradius 

no 

r g 

normalized gyroradius 

S 

length of particle orbit in 
tersecting the aperture 

T z 

axial temperature 

T 1 

transverse temperature 

V 

retarding potential 

v x 

component of velocity 


along x-axis 
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V 

y 

component of velocity 
along y-axis 

v z 

component of velocity 
along z-axis 

v z 

normalized axial velocity 

V 1 

transverse velocity 

V I 

normalized transverse 
velocity 

r 

constant, defined in 
eq. (35) 

6 

transverse distance 
particle travels while 
traversing the aperture 
length (B = 0) 

V 

transmission function 

V 

normalized transmission 
function 

% 

transmission function 
(B = 0) 

e 

orbital angle 

k 

normalized distance 


p distance from center of 

aperture to guiding 
center of particle 

p normalized distance from 

center of aperture to 
guiding center of 
particle 

p^ lower limit, eq. (26c) 

Pq distance from center of 

aperture to guiding 
center of particle 
(B = 0) 

r 

p lu upper limit, eq. (26a) 

P 2u upper limit, eq. (26b) 

cr surface density of guid- 

ing centers 

r time needed for particle 

to traverse the aper- 
ture length 

(p angle defined in fig. 3 
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APPENDIX B 


LIMITS OF INTEGRATION 
Upper and Lower Limits for Class II Particles 

Figure 9 shows the upper and lower limits on p for a typical class n particle at 
various values of 0. As previously mentioned, particles with P > P u or P < P^ will 
make the integrand (S - r^d) negative, and are thus absorbed by the aperture wall. Refer- 
ring to figure 10 and applying the law of cosines gives 



Figure 9. - Upper and lower limits for typical class II particle as function of orbital angle. 



(a) Upper limit. (b) Lower limit. 

Figure 10. - Upper and lower limits on distance from center of aperture 
to guiding center of particle for class II particles. 
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Solving equation (Bl) for p u ^ and using the notation that has been used in this report 
result in 


p 2u = r g cos (f) + y <R 2 - r g si " 2 (f) < B2a > 

- % c ° s (i) - V^RS <B2b) 

The upper and lower limits coalesce when 

R = r g sin ^ (B3) 

This implies that the integral vanishes for values of 6 or r such that 

B 

r g sin (i) > r < B4 > 


Upper Limit of Integration for Class lb Particles 

The upper limit of integration on p for typical class lb particles at various values 
of 6 are shown in figure 11. The upper limit is obtained in the same manner as the 
limits were obtained for class II particles and is given by 

Pm - r g POP (f) + |/r 2 - r| sin 2 (B5) 



Figure 11. - Upper limit on distance from center of aperture to guiding center 
of particle p for typical class lb particle at various values of orbital angle 0. 
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APPENDIX C 


INTEGRATION OF THE TRANSMISSION FUNCTION 
Integration of Equation (12) 


It was previously noted that the problem of finding a closed form expression for the 
transmission function consists in evaluating an integral of the form 


1 = 


f [» -.-I 


P dp 


J* 2 cos 1 |ajp - — ^ p dp - 


(Cl) 


Consider the equivalent expression for the integral in equation (Cl): 


= / b 


-1 


V = / x cos u dx 
a 


(C2) 


Integrating by parts converts this equation to 

lb „b 


I’ = ^ cos- 1 u 


+ 

a 'a 


/ 




du 

dx 


dx 


(C3) 


Hence, the original integral becomes 
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I HIM ■ ■ 


Concerning ourselves with the remaining integral 


where 


f a i + | P 2 d P 


2 1 1 / 2 


/'* tv .-, 

J |" a i p + (1 + 2a^ 


)p dp 


Q v 2 2 

a 2 )p -a 2 


nl/2 


i/ 


(ajV + a 2 ) dv 


= ^1 f h vdv + ^2 f 

2 J 2 Vv 2 V2 VV 

d ct 


dv 


2a, 


ib‘ 


+ 1 ; 4a l a 2 sln -l l 2 *‘l v ~ 1 - 2a l a 2l 


4a i 

2 1 


yi + 4a 


(C5) 


V = 


v= P 


-a^v 2 + (1 + 2aja 2 )v - a 2 


(C6a) 

(C6b) 


And finally 

r 

j , 2 

I =< p cos 




2 |^ -a iP + (1 + 4a ^ a 2 )p - a 2 J 


1/2 


2a, 


2 2 


1 + 4ai&2 cm" 1 -1 ~ 2 ^ 2 ' 


4a, 


^/l + 4a^£ 


(C7) 


28 


Integration of Equations (5) and (6) for the Case e = 0 


In deriving the transmission function, the statement was made that the number of 

o 

particles incident on the aperture was represented by onR , independent of the value of 
the gyroradius. The proof follows. 

Taking the case with r g less than R, equation (5) may be written as 


N thru = an( - R ~ r g> + CTl 

Substituting the value of I given in equation (C7) 


R+r 


g 


R-r 


(C8) 


g 


N thru = CT7T ( R " r g> + 


p 2 cos 




a 2\ + (1 + 2a i a 2^ - a£^ 


1/2 


2a, 


2 2 

1 + 4a^2 j 2a^p 1 — 


4a‘ 


^1 + 4aj£ 


J 


R+r 


(C9) 


R-r 


g 


Substituting the values of a^ and ag into equation (C9) one obtains 


N thru = CT ^ R - V +a < 


,2 2 d 2 > 

2 -1 ( p + r g ' R 

p cos 1 § - 


2pr 


g 


= ott(R - r g ) + ct 


= aTi(R - r g ) + a 


+ R sin 


(R + r g ) 2 cos 1 (1) - (R - r g ) 2 cos 1 (-1) 


4 2 o 2 (R 2 - r 2 ) 

-P_ + (R 2 + r 2 )e_-V SL 

4 g 2 4 


2 2 2N 

2 . -i IP ~ r ~ -R 


1/2 


_S_ 

2r R 
g 


iR+r 


R-r 


- 0 + 0 + R 2 sin -1 (l) - R 2 sin _1 (-l) 


x2 ttR 2 77 R 2 
-^(R - r ) + + 

g 2 2 


= ottR 


(CIO) 
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For the case where r is greater than R, equation (6) becomes 


|H +R 

N — or I & 

" to 1 r -R 
g 


= a|(R + r g ) 2 cos" 1 (1) - (r g - R) 2 cos -1 (1) - 0 + 0 + R 2 sin" 1 (1) - R 2 sin" 1 ( -l)j 


= ctttR 


Thus for both ranges of r , the total number of particles impinging on the front sur- 

& 2 

face of the aperture comes out the expected value of oxR even when expressed in this 
complex form. 
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APPENDIX D 


COMPUTER PROGRAMS 

The programs that follow are written in standard FORTRAN IV language and do not 
require any special functions other than those normally supplied. 


Zero Magnetic Field 

A flow diagram for the program which calculates the distribution function for the 
zero magnetic field case (APERTO) is given in figure 12. This program varies the 
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rm 


normalized axial energy (denoted by QV within the program) from 0 to 5 in steps of 0. 1 
and computes the aperture exit distribution function (eq. (34)) for each value of QV. 

The integration over the transverse distribution function is performed by subroutine 
SIMSN, which is a straightforward Simpson’s 1/3 rule numerical integration. Although 
the limits are supposed to vary from 0 to «>, the upper limit is arbitrarily restricted 
to a maximum value of 5 (exp(-25) is after all a rather small number). 

A flow diagram for subroutine SIMSN is given in figure 13. This subroutine is de- 
signed to repeatedly double the number of steps used in the integration until two succes- 
sive values of the integral agree to within a specified number of digits. 



Figure 13. - Flow diagram - 
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The complete program listings follow, 


COMPUTATION OF THE DISTRIBUTION FUNCTION FOR THE EFFECT OF AN 
APERTURE ON A TWO TEMPERATURE MAXWELLIAN OISTRIBUTION WHEN B=Q. 


KK = 0 

PRINT OUT HEADINGS* 

10 WK I TE i 6. 6 10 ) 

610 FORMAT! iHl . 10X .28HEFFEC T OF APERTURE WHEN 8*0.///) 

WRI TEI6, 6001 

) FORMAT! 1H ♦ 10X .2H0V.6X. 9HD i ST . FN* .6X.3HL/R.6X. 5HTP/TZ/// ) 

INITIALIZE PROGRAM AND SUPPLY REQUIRED CUNSTANTS. 

COMMON GAM A, PI 
P 1=3* 141 593 

SUPPLY INPUT DATA. 

XLR IS THE APERTURE LENGTH TO RADIUS RATIO. 

TPZ IS THE TRANSVERSE TO AXIAL TEMPERATURE RATIO. 

READ! 5. 50 0) XLR . TPZ 
500 FORMA T ! 2F 10 .3 ) 

C 

C END OF INPUT DATA. START OF CALCULATION OF DISTRIBUTION 
C FUNCTION. 

C 

DIMENSION DISTFN! 51J.0VV! 511 
DO 200 J =1. 51 
XJ=J- I 
QV=XJ/10. 

OW! J 1 =0 V 
VZ=SURTi QV1 

GAMA=2.828428*VZ/1 SGR T ! TPZ ) *XLR) 

IF ! GAMA .ED. O.OJ GO TO 210 
GO Til 220 

210 DISTFN! J J *0.0 

GO TO 200 
220 A=6XP!-QV) 

IF i A .LT. 1.0 E- 30) A=0. 0 

EXTERNAL ARG 

DlFF=0.0 

IF ( GAMA .LE. 5.) b = SI MSN ( 0. 0 . GAM A ,4 . ARG. 01 FF 1 
IF ( GAMA .GT. 5.) B*SI MSN ! 0. 0 . 5. .4 # ARG • D IFF ) 

IF ! B .LT. 1.0 E-30) B=0.0 
DISTFNIJ) =A*B 
200 CONTINUE 

DO 300 N= l * 51 

WRITE! 6, 605 10 WIN I .01 STFN(N) .XLR. TPZ 
605 FORMAT! 1H . F13 . 2. E 1 5. 4, 2F 1 0. 2 1 
300 CONTINUE 
KK=KK + 1 

IF IKK .LT. 99) GO TU 10 

STOP 

END 


M 10 
H 20 
H 30 
M 50 
M 60 
M 70 
M BO 
M 90 
M 100 
M 102 
M 104 
M 110 
M 120 
M 130 
M 140 
M 150 
M 160 
M 170 
H 180 
M 190 
M 200 
M 210 
M 220 
M 230 
M 240 
M 250 
M 260 
M 270 
M 280 
M 290 
H 310 
M 320 
M 330 
M 340 
M 350 
M 360 
M 365 
M 366 
M 367 
M 368 
M 369 
M 370 
M 380 
H 390 
M 395 
M 400 
M 40 5 
M 410 
M 420 
M 430 
M 470 
M 580 
M 590 
M 610 
M 615 
M 620 
M 630 
M 640 


SIBFTC SUB1 


c A 10 

C THIS FUNCTION IS THE TRANSVERSE DISTRIBUTION FUNCTION TIMES ETA. A 20 

c A 30 

FUNCTION ARG! Y ) A 40 

COMMON GAMA. PI A 50 

X=2.* Y/GAMA A 55 

ETA=I2.*ARC0S! X/2. > -X*SORT I 1 . 25 *X*X> )/PI A 60 

ARG=Y*ETA*EXPI-Y*Y/2.) A 70 

RETURN A 80 

END A 90 



Ill IHII 


C 1 

C FUNCTION SIMSN I A.B ,N * FN , D IFF 1 2 

C 3 

C SIMPSON 1/3 RULE INTEGRATION* 4 

C A IS THE LOWER LIMIT* B THE UPPER LIMIT • N GIVES THE NUMBER OF 5 

C SIGNIFICANT FIGURES THE ANSWER CAN BE EXPECTEO TO BE CORRECT TO, 6 

C AND 100*0 IF F IS THE PERCENT DIFFERENCE BETWEEN THE LAST TWO 7 

C COMPUTED VALUES OF THE INTEGRAL. 8 

C 9 

C THE CALL SEQUENCE IS 10 

C 11 

C EXTERNAL FN 12 

C DIFF=0.0 13 

C ANS = SIMSNIA*B*N*FN#D1FF ) 14 

C 15 

FUNCTION SIMSNIA.fi, N V FN,DIFF) S 10 

EXTERNAL FN S 20 

K=0 S 30 

8N=10* S 40 

ERK0R*1. S 50 

DO 100 J = 1 , N S 60 

100 ERROR -ERROR /1 0. S 70 

200 H= ( B— A 1 / HN S 80 

EVEN^O.O S 90 

ODD=0*0 S 100 

NA=HN-l. S 110 

NB = H\I— 2. S 120 

DO 300 KA-1.NA,2 S 130 

XKA-KA S 140 

X=A+XKA*H S 150 

300 000=OQD*FNIX) S 160 

DO 400 KB=2,NB,2 S 170 

XKB=KB S 180 

X*A+XK8*H S 190 

400 EVEN=E VEN+FN l X ) S 200 

IF ( K .GT* 01 GO TO 500 S 210 

K-K + 1 S 220 

ANSi=H*(FN(A)+FN<B)+4. *000+2. *EVEN) /3. S 230 

HN = 2-*HN S 240 

GO TO 200 S 250 

500 IF(K .GT, 11 GO TO 600 S 260 

K-K+l S 270 

ANS2=H*(FN(AI«-FNlB) + 4.*0DD+2. *EVEN) /3- S 280 

IFIANS2 *E0. 0.01 GO TO 700 S285 

0IFF=A8SI ( ANS 1-ANS2 1 /ANS2 i S 290 

IFIDIFF *LT. ERROR 1 GO TO 700 S 300 

550 H N=2.*HN S 310 

IF ( HN *LT. 100*1 GO TO 200 S 320 

600 ANS1=ANS2 S 330 

ANS2=H*I FNIAI +FN( B 1 *4.*ODO+2* *E VEN) /3- S 340 

DIFF=ABS( I ANS1-ANS21/ANS2) S 350 

IFIDJFF .IT. ERROR) GO TO 700 S 360 

GO T3 550 S 370 

700 SIMSN*ANS2 S 380 

RETURN S 390 

*=*0 S 400 


Nonzero Magnetic Field 


The equations used in evaluating the distribution function for the case of a finite 
magnetic field are (25b) to (25e), (26) and (27). The flow diagram for the main program 
(APERT) is given in figure 14. The calculation of the distribution function proceeds via 
subroutines FNI and FNA. 

Subroutine FNI solves the relation 


FNI = (1 



e 


-V 


2 

1 


/ 2 


r>~> 



(Dl) 
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Figure 14 - Flow diagram - program APERT. 
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and is appropriate when the value of v z is such that only class la particles can pass 
through the aperture. 

Subroutine FNA calculates eta times the transverse distribution function for the 
general case and is depicted in figure 15 . 



Figure 15. - Flow diagram subroutine FNA. 
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Figure 16. - Flow diagram subroutine XL 


Figure 16 is a flow diagram for subroutine XI, which is a solution of the integral 
used in obtaining eta (eq. (25e)). The program listings follow. 
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C ^ 10 

C COMPUTATION OF THE EFFECT OF AN APERTURE ON A TWO TEMPERATURE M 20 

C MAX WFLl I AN DISTRIBUTION FUNCTION. M 30 

C M 40 

C A MAGNFTIC FIELD IS ASSUMED PERPENDICULAR TO THE APERTURE PLANE. M 50 

C M 60 

C THF NORMALIZED APERTURE LE NG TH=4 7200*B*L/ SORT ( TZ I M 70 

C B IS THF MAGNETIC FIELD IN TESLA. M 80 

C I IS THF APERTURE LENGTH IN METERS. M 90 

C T? IS THF AXIAL TEMPERATURE IN E.V. M 100 

C M 110 

C PRINT OJT HEADINGS AND SET COUNTER. M 120 

C M 130 

KK - 0 M 140 

10 HRITFIo.610) M 142 

610 FORMAT ( lHl .20X.43HEFEECT UF APERTURE UN DISTRIBUTION FUNCTION///) M 144 

WRITF(6.bOO> M 150 

600 0 FnRMA T ( 1 H • 5X . 2H J V . 6X . 9HD 1ST. F N. ,6 X . 3H L/R . 1 1 X , 2 HLN * 1 2 X . 2 H RN , M 160 

l I2X.5HTP/T Z///J M 16L 

C M 170 

C INITIALIZE PROGRAM AND SUPPLY RE UU I R ED CONSTANTS. M 180 

C M 190 

RFAL LN* LR M 200 

COMMON L N « Lk * TZTP.PI , RN * VZ M 210 

PI = 3. 141593 M 220 

KK =K< + 1 M 230 

C M 240 

C INPUT DATA FOLLOWS. M 250 

C LR IS THE LENGTH TO RADIUS -RATIO UF THE APERTURE. M 260 

C LN IS THE NORMALIZED APERTURE LENGTH. M 270 

C T Z TP IS THE RATIO UF THE AXIAL AND TRANSVERSE TEMPERATURES. M 280 

C M 290 

RFAD1 5.500) LR.LN.TZTP M 300 

50D FORMAT I 3F10.4) H 310 

TPTZ=l./T7TP M 320 

C M 330 

C END OF INPUT DATA. START OF CALCULATION OF DISTRIBUTION M 340 

C FUNCTION- M 350 

C KN rs THh NORMALIZED APERTURE RADIUS. M 360 

C M 370 

KN=2.*PI *LN*SQRT(2.*TZTP)/LR M 380 

DIMENSION D IS TFN I51J.QVVI51) M 390 

DO 200 J= 1 , 51 M A 00 

XJ=J- 1 M 410 

QV=XU / 10. M 420 

OVVIJJ-OV M 430 

VZ = S JR T ( UV ) M ^40 

C M <*50 

C CHECK T) SEE IS VI IS GT THAN MIN VALUE RORD. M 460 

C M A 7 0 

I F ( V/ .GT- LN) GJ TO 100 M 480 

t RN-RN M 490 

EXTFRNAL FN l M 500 

DIFF^O.O M 505 

IFIRN -LT. 5.) DISTFNI J)=fcXPI-UV> *SI MSN ( J . 0 , RN, 4 , FNI . DI FF > M 510 

IF ( R M -GE. 5.) DISTFMJ)-EXP( — OVI *SIMSN(0.0»5-.4,FNI»DIFF) M 515 

GO TD 200 M 520 

100 CONTINUE M 530 

FXTERNAL FN A M 5^0 

01 FF = 0 .0 M 545 

DISTFNU ) =EXP I-JV)*SI MSNI0.G.5. .5 .ENA.OIFF) M 550 

200 CONTINUE M 560 

DO 3D0 N-L . 51 M 600 

WR I TF I 6* 60 5 I OVVIN).DISTFN(N) .LR.LN.RN.TPTZ M 710 

605 FORMA TI 1 H . F8 . 4 . 5E 1 4. 4 ) M 720 

300 CONTINUE M 730 

IF ( K< -LT. 175) GO TO 10 M 740 

STOP M 750 

END M 760 
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A 10 

THIS FUNCTIflN IS THE TRANSVERSE DISTRIBUTION FUNCTION TIMES A 20 

FT A FOR CLASS IA PARTICLES. A 30 

A 40 

FUNCTION FNI(VT) A 50 

COMM IN LN.LR, T2TP.PI ,RN,VZ A 60 

REAL LN.LR A 70 

IFIVT .GT. 7.0) GO TO 10 A 72 

GO TO 20 A 73 

10 FN 1=0.0 A 74 

GO TO 100 A 75 

20 RG=VT/RN A 80 

FN I =v/ T« 4 l.-RG) *( l.-RG)*EXP I-VT*VT/2. ) A 90 

100 CONTINUE A 95 

RETURN A 100 

FND A 110 


1 

FUNCTION S I MSN ( A * 8 • N • FN , 0 I FF ) 2 

3 

SIMPSON 1/3 RULE INTEGRATION. 4 

A IS THE LOWER LIMIT, ti THE UPPER LIMIT , N GIVES THE NJMBER OF 5 

SIGNIFICANT FIGURES THE ANSWER CAN BE EXPECTED TO BE CORRECT TO, 6 

AND 1 00*0 IFF IS THE PERCENT DIFFERENCE BETWEEN THE LAST TWO 7 

COMPUTED VALUES OF THE INTEGRAL. 8 

9 

THE CALL SEOUENCF IS 10 

11 

FXTERNAL F.M 12 

DIFF=0.0 13 

ANS = SIMSN( A .B.N.FN.DIFF ) 14 

15 

FUNCTION SIMSN(A,B,N»FN,D1FF) S 10 

EXTERNAL FN S 20 

K=0 S 30 

HN = 10. S 40 

ERROR = 1. S 50 

DU 100 J = 1 • N S 60 

100 ERROR = FRRDR/IO. S 70 

200 H= ( d- A ) / H.N S 80 

FV FN = 0.0 S 90 

nOD=(). 0 S 100 

NA = HN- 1 . S 110 

N B = HN —2m S 120 

DO -10 0 <A=1,.NA.2 S 130 

XK A A S 140 

X=A + XKA*H S 150 

300 UDD=ODO+FN( X) S 160 

DO 400 <H=2.NB,2 S 170 

XKB=KR S 180 

X= A + XKB* H S 190 

400 EVEN = EVEN+FNI X ) S 200 

I F ( K .GT. 0) GO TO 500 S 210 

K=K+1 S 220 

ANS 1 = H*( FN( A) *-FNl B ) + 4. +0DD + 2 • #EVEN) / 3. S 230 

HN = 2 . * HN S 240 

GO TO 200 S 250 

50 u IF ( K .GT. 1) GO TD 600 S 260 

K=K+l S 270 

ANS2=H*< FN( A) +FNI B)+4. *000+2. *tVEN)/3. S 280 

DIFF=ArtS( ( ANS1-ANS2J/ANS2) S 290 

IF ( D I FF .LT. ERROR) GO TO 700 S 300 

550 HN = 2. * HN S 310 

I F ( HN .LT. 100.) GO TO 200 S 320 

.600 AN SI = AN S? S 330 

ANS2=H*( FN( A) +FN( B )+4.*0DD+2.*EVEN)/3. S 340 

DIFF-ABSI t ANS1-ANS2) /ANS2) S 350 

IFIDIFF .Lf. ERROR) GO TO 700 S 360 

GO TO 550 S 370 

70 0 S I MSN =AN S2 S 3 80 

RETURN S 390 

FND S 400 
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B 10 

THIS FUNCTIUN IS THE TRANSVERSE DISTRIBUTION FUNCTION TIMES B 20 

ETA FOR ALL PARTICLES WHOSE AXIAL VELOCITY EXCEEDS THE MINIMJM B 30 

REQUIRED. B 40 

B 50 

RfllJ IS THE UPPER LIMIT OF INTEGRATION. B 60 

ROL IS THE LOWER LIMIT OF INTEGRATION WHEN RG .GE. 1. B 70 

B 90 

FUNCTION FN At V T ) B 100 

COMMON LN.LR. TZTP.PI ,RN,VZ B 110 

REAL L N * L R B 120 

EXTERNAL XI B 130 

RG=VT/KN B 140 

JFtRG .EQ. 0.0) 30 TO 35 B 141 

Gfl H 40 B 142 

35 F TA= l - B 143 

GO T) 100 B 144 

40 THFTA=2.*P I+LN/VZ B 160 

IF I THETA .GT. 2.*PI) THETA=2.*PI B 170 

25 RQA=RG*COS< THETA/2. ) B 200 

1FIRG .LT. 1.) GJ TO 30 B 202 

I F ( THE TA .GT. ? .* AR SI Ml . /RG ) ) GO TO 10 B 204 

30 ROB=L.-RG*RG*SIN( THETA/2. )*SI NITHETA/2. ) 8 210 

IF ( ROB .LE. 0.0) GO TO 10 B 220 

GO TO 20 B 230 

10 FT A =0 . 0 B 240 

GO TO 100 B 250 

20 RORA-SQR T ( ROB ) B 260 

RO U= R O A + R 0-B A B 270 

ROL =R0 A-ROB A B 280 

0 IFIR3 .LT. l.| ETA = XI ( VT , RGU* RN)/PI+.5 B 300 

1 +THETA*t l.-RG) *t l.-RG) /(2.+PI) B 301 

IF ( RG .GE. 1.) ETA=(XIIVT,RQU* RN) -XI I VT »RQL»RN) ) / P r B 310 

I F ( ETA .LT. 0.0) ETA = 0.0 B 320 

IF ( ETA .GT. 1.) E TA = 1. B 325 

100 FNA=VT*ETA*EXP(-VT*VT/2. ) B 340 

RF TURN B 350 

END B 360 



FUNCTION Xi<VT.RO,RN) 




10 

C 





20 

c 

THIS FUNCTION GIVES THE VALUE OF THE INTEGRAL USED TO CALCJ LATE ETA 


30 

c 

VT IS THE NORMAL l 7 E 3 TRANSVERSE VELOCITY, 

RO THE 

LIMIT OF 


40 

c 

INTEGRATION. TFETA THE URBITAL ANGLE, AND 

RN THE 

NORMALIZED 


50 

c 

APERTURE RADIUS. 




60 

c 





70 


CnMMlN LN.LR, TZTP.PI ,RN,VZ 




80 


REAL LN.LR 




90 


RG=VT/RN 




100 


IF ( RG -EJ. 0.0) GO TO 50 




1 10 


RGS=RG*RG 




120 


A=( R7+R.7-KGS- 1 . 1 / f 2.*RG) 




130 


IF ( A .LT. 0.0) XI A=-AKSIN(-A) 




140 


I F ( A -GE. 1.) XIA=PI/2. 




150 


I F ( A .GE. 0.0) XI A=AR SI N ( A ) 




160 

70 

B=(R1*R0+RGS-1. ) / (2.*R0*RG> 




230 


IF ( B .GE. 1.) X I B = 0 • 0 




260 


IFI-B -GE. 1.) XI B =RD*RD*P 1 




270 


IFIABS(R) .LT. 1.) X I B =RO*RO* ARC OS t B ) 




280 

45 

C=-RJ**4/4.+R0*R3*I R3S+ 1. ) /2 - - t 1 . -RGS) ** 2 /4. 



340 


I F ( C .LT. 0.0) C = 0.0 




355 


XIC=SQRT I C ) 




360 


THFTA=2.*PI *LN/V/_ 




380 


IF t THFTA .G1. 2.* PI) THETA=2.*Pl 




390 

16 

XI =X I A+X IB- XI C — THE TA*R0*R0/2. 




420 


GO TO 60 




430 

50 

XI=0.0 




440 

60 

CONTINUE 




450 


RETURN 




460 


END 




470 
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